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ABSTRACT 

The connection between the two-point and the three-point correlation func- 
tions in the non-hnear gravitational clustering regime is studied. Under a scaling 
hypothesis, we find that the three-point correlation function, obeys the scal- 
ing law C 0^ ^ 2m+2™ in the nonlinear regime, where ^, m, w, and e are the 
two-point correlation function, the power index of the power spectrum in the 
nonlinear regime, the number of spatial dimensions, and the power index of the 
phase correlations, respectively. The new formula reveals the origin of the power 
index of the three-point correlation function. We also obtain the theoretical 
condition for which the "hierarchical form" ^ oc is reproduced. 

Subject headings: cosmology: theory - large-scale structure of universe - meth- 
ods: - analytical - methods: - numerical - methods: statistical 

1. Introduction 

The formation of large-scale structures in the Universe is one of the most important 
and interesting problems in cosmology. It is generally believed that these structures have 
developed by a process of gravitational instability from small initial fluctuations in the density 
of a largely homogeneous early Universe. Hence it is very important to clarify the physical 
mechanism of the evolution of density fluctuations. Here we consider the density fluctuations 
of collisionless particles, such as dark matter, because our interest is concentrated on the 
effects of self-gravity. 

The inflationary scenario predicts a random Gaussian fleld possessing the properties of 
statistical homogeneity and isotropy as the primordial fluctuations (Guth 1981; Guth & Pi 
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1982; Bardeen et.al. 1986). The statistical properties of the Gaussian field are completely 
characterized by the two-point correlation function (hereafter 2PCF) or the power spectrum, 
since higher order correlations vanish, while only the 2PCF has a certain value (Adler 1981). 
However, non-linear gravitational clustering gives rise to non-zero values of the higher-order 
correlations, even if the primordial density fluctuations were set Gaussian. This is reasonable, 
since phases between Fourier modes in the bispectrum are considered to be distributed non- 
randomly but correlated in some ways in the non- linear regime, where galaxies or galaxy 
groups are strongly clustered around each other (Scherrer, Melott & Shandarin 1991; Goles 
& Ghiang 2000; Ghiang 2001; Ghiang, Naselsky & Goles 2002; Watts & Goles 2003; Watts, 
Goles & Melott 2003; Matsubara 2003). Hence the higher order statistics, hke the three- 
point correlation function (hereafter 3PGF). become essential when the evolution of the 
density flelds is governed by non-linear gravitational clustering. So it is quite important for 
understanding non-linear evolution to investigate the properties of the 3PGF. 

Generally, it is hard to deal with the exact formula of the 3PGF both theoretically and 
numerically, since it is much more lengthy and complicated than 2PGF. Historically, however, 
Peebles & Groth (1975); Groth & Peebles (1977) had proposed a simple assumption: the 
3PGF can be expressed by the product of the 2PGF as 

Cabc — Qi^ab^bc + ^bc^ca + C ca^ab ], (1) 

where ^, and Q are the 3PGF, the 2PGF, and a certain constant, respectively. The 
expression (1) or ( oc ^'^ is often called the "hierarchical form" (Pry 1984). The data from 

the Zwicky and Shane- Wirtanen catalogs arc in good agreement with this assumption when 
we choose Q ~ 0.85 for Zwicky and Q ~ 1.24 for Shane- Wirtanen (Peebles & Groth 1975; 
Groth & Peebles 1977). 

However, there is a serious fault with the proposal (1). It has no theoretical ground, 
although some observational or numerical data may be wcll-fltted. One explanation may be 
that Eq. (1) for ( with Q = 1 coincides with the Kirkwood superposition approximation 
which is familiar in liquid physics and turbulence theory, if ^ah^bc^ca ^ 1 (Davies & Peebles 
1977; Ichimaru 1973). However, this approximation is not appropriate for the non- linear 
regime, since the distribution of the gravitational sources is strongly correlated, ^ S> 1, in 
such regimes. Hence this resemblance is helpless to compensate for the theoretical defect. 

Indeed, there were many theoretical investigations of the 3PGF (Fry 1984; Ruamsuwan 
& Fry 1992). Nevertheless most of the analyses have been based on the "hierarchical form" 
assumption. The coefficient Q under the "hierarchical form" was estimated by assuming self- 
similarity (Pry 1984). The BBGKY hierarchy was analyzed by assuming the "hierarchical 
form" as well as self-similarity (Ruamsuwan & Pry 1992). In the first place, however, there is 
no reason why the 3PGF should necessarily depend only on the second power of the 2PGF, 
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even if we accept the self-similarity. Indeed, Yano & Gouda (1997) suggests the possible 
existence of solutions that do not satisfy the "hierarchical form" on the basis of the BBGKY 
equations. Therefore an investigation of the physical adequacy of the "hierarchical form" 
itself and more reliable and useful formulae for the 3PCF are pressingly required. 

In this paper we study the connection between the 2PCF and the 3PCF in the non- 
linear gravitational clustering regime. We analyze these functions using a scaling hypothesis. 
It is expected that the time evolution of the statistics obeys also some self-similar rules 
fundamentally, since gravity is scale free (Davies & Peebles 1977; Fry 1984; Ruamsuwan & 
Fry 1992; Yano & Gouda 1998). Many authors have investigated the self-similarity of the 
2PCF or the power spectrum. Recent works using N-body simulations (Colombi, Bouchet, 
& Hernquist 1996; Couchman & Peebles 1997; Jain, Mo & White 1995) have shown that 
self-similar evolution of the power spectrum can be satisfied when the initial power spectrum 
is scale free, and the power spectrum in the non-linear regime obeys some scaling laws. The 
scaling hypothesis itself is very familiar. However, the unique proposition in this paper is 
that we apply the scaling hypothesis to the bispectrum as well as to the power spectrum, 
and we never assume the "hierarchical form". Consequently, we succeed in deriving a new 
formula for the 3PCF, focusing on the power indices. The relevance to the "hierarchical 
form" is also discussed. 

Our paper is organized as follows. In Sec. 2 we introduce the general formula of the 2PCF 
and the 3PCF. We obtain in Sec. 3 a new relationship between these functions the under the 
scaling hypothesis. In Sec. 4 we calculate the 2PCF and the 3PCF in the non-linear regime in 
the one-dimensional Einstein-de Sitter model numerically, in order to confirm the analytical 
results. The final section is devoted to summary. 



In this section we consider the 2PCF and the 3PCF. Let us denote by p the mean density 
and take p{x) to be the density at a point specified by the position vector x with respect to 
some arbitrary origin. As usual, we define the fiuctuation to be 
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General formula for the correlation functions 



5{x) 



p{x) - p 



(2) 



P 



We assume this to be expressible as a Fourier series: 




(3) 
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and 

'^(^^ = (If / ^k^"^''"d'"k, (4) 

where 6j^, w and L are a Fourier amplitude, the number of the spatial dimensions, and the 
length of the side of the box, respectively. The 2PCF is defined in terms of the density 
fluctuation by 

= (( / [ SA'--'**"'' ■ a^e-*-'-rf»H»i') J , (5) 

where the angular brackets in Eq.(5) are expectation values, formally denoting an average 
over the probability distribution of 5{x). The first average, ( )x, means simply a volume 

normalization over a sufficiently large patch of the Universe. The second average, ( ), is 
interpreted as that across independent samples in many Universes. Applying this rule and 
performing the volume integration gives 

^ J^I P(^>"^'"d'"k, (6) 

where 5c is the Dirac delta-function and P{k) is the power spectrum. Similarly, it is possible 
to calculate the 3PCF, which is defined as 

C(r,s) = {{6{x)6{x + r)6{x + s))x) 

The volume integration gives 

^(^'*) = JJJ{^k^k'^k'')^D{k + k' + k'')e-'^^-''+'^'-^^rM 



L 



2w 



(27r)2«' 



{^kV-k-k')^~'^^'''^^''^^dVk' 
J j B{k, k')e-'^^ '^+^'-^U'"k(Fk', (8) 



where B{k,k') is the bispcctrum (Fry 1984). The vohimc integration gives rise to the 
triangular constraint by the Dirac delta-function, which is similar to Eq.(6) (see also Fry 
(1984); Watts & Coles (2003); Matsubara (2003)). 
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3. The connection between the 2PCF and the 3PCF under the scaUng 

hypothesis 

In this section wc derive the connection between the 2PCF and the 3PCF, under the 
scahng hypothesis. This assumption is adequate, since we consider that fluctuations evolve 
purely due to gravity. First, we assume that the power spectrum in the non-linear regime 
obeys a power law as 

p{k) = (i^r) 

oc yfe"*, (9) 

where the index is a value in the non- linear regime. Substituting Eq.(9) into Eq.(6) gives 

e(r) oc j ¥^(rk 

oc r-^^+^\ (10) 
Next, the bispectrum in Eq.(8) is decomposed by 

B{k,k') = (|5fe||5fe,||5_fe_,,|e^(^'=+<^^'+<^--')> , (11) 

where (f)k is the Fourier phase of the modes with wave number k. Now we consider the Fourier 
phases in Eq.(ll) in the non-linear regime. Let us consider that the initial fluctuations are a 
random Gaussian field whose phases arc distributed randomly. Induced by the gravitational 
instability, however, the density fluctuations become stronger and stronger, which means 
the phases between Fourier modes become distributed non-randomly (Scherrer, Melott & 
Shandarin 1991; Coles & Chiang 2000; Chiang 2001; Chiang, Nasclsky & Coles 2002; Watts 
& Coles 2003; Watts, Coles & Melott 2003). That is, the Fourier phases are also correlated 
due to the gravitational clustering. It is not difficult to imagine that the phases synchronize 
in the case of the gravitational collapse, like a pancake, for example. 

Generally, the 3PCF or the bispectrum depends on the shape of the configuration. If we 
fix the shape of the triangle, however, the 3PCF and the bispectrum are determined by one 
parameter. Here we assume the bispectrum in the strong non-linear regime obeys a scaling 
law as well as the power spectrum. Hence the expectation value of the Fourier phase in the 
bispectrum (11) is described as 

l^Qi{4>kt+4>k't+<t> -kt-k't)'^ oc t^'^ {e^^'t'i=+^ki+4>-k-ki)^ (^\2) 



in the strongly non-linear regime. The scaling hypothesis for the Fourier phases corresponds 
to the following two characteristic possibilities physically: first, the Fourier phases remain 
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correlated, and the expectation value converges to some non-zero constant value through 
non-linear gravitational clustering. In this case the value of a in Eq.(12) is zero. Second, 
the Fourier phase correlations are finally relaxed due to non-linear gravitational clustering, 
although the phases may be correlated once by gravitational instability. In this case the 
expectation value of the Fourier-phase correlations converges to zero, and the value of a in 
Eq.(12) is some finite value that depends on the degree of convergence. That is, the constant 
a should be a non-negative value. Similarly, other terms in Eq.(ll) are also assumed to obey 
some scahng laws (For a detailed description, see Appendix A). 

From Appendix A, the 3PCF is reduced to 
C(r,s) oc jj B{k,k')d''kd'^k' 

+ •••. (13) 

Where j3 and 7 are the scaling factors of variance and other constants. It is noting that each 
value of these factors should be non-negative so that the 3PCF should not diverge. Then 
the 3PCF can be approximated by 

C(r,s) oc ^"^^"^^'""^^C , (14) 

C(r) oc r-^lrn+2w-e) ^ (^5) 



that is, 
where e is defined as 



e = mm(Q;, 72y, /3 73y, /3 + 742/)- (16) 
From Eqs.(lO) and (15), we obtain a new expression to connect between C, and ^ as 

3m+4w-2e 

( OC ^ 2m+2t<, . (17) 

The proportional expression (17) states that the 3PCF is proportional to the power of the 
2PCF. The power index depends on the power index of the power spectrum in the non-linear 
regime, m, and the scaling factors of the bispectrum, e. The formula clarifies the source of 
the power index of the 3PCF under the scaling hypothesis. 

Now we consider the case a 7^ 0, where the phase the phase correlations get loose by the 
non-linear gravitational clustering, although they are tied tight once. Still the information 
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of the convergence survives in the power index of the formula (17) . We can know the power 
index of the 3PCF by calculating the power index of the power spectrum and the phase 
correlation in the bispectrum (11). In particular, if the equality e = —m/2 is satisfied, the 
"hierarchical form", C, oc can be reproduced. For example, we consider the case ^ oc r~^-^ 
[Davies&Pccblcs 1977]. In this case, the index m is equal to —1.2, and our formula agrees 
to the "hierarchical form", if e = 0.6. 

Next, we consider the case where the Fourier phases keep correlated through the non- 
linear effects. In this case the expectation value of the phase correlations converges to a 
finite constant value, « = 0. Then the power index of the formula (17) depends only on the 
index of the power spectrum and the number of the spatial dimension. 



Accordingly, the 3PCF can be described as the closed form only by the 2PCF. As a typical 
example, it is expected that the relation is realized at the moment of the gravitational collapse 
into a pancake, since the phases are synchronized strongly. In addition, if the power index 
of the power spectrum is zero, m = 0, the 3PCF is proportional to the second power of the 
2PCF. Recall that the index of the power spectrum is around zero in the transitional region 
from the linear regime to the non-linear regime, if the power index of the power spectrum of 
the initial random Gaussian field is not negative value. Considering these facts, we suppose 
the "hierarchical form" might be also realized accidentally around such quasi-linear scales. 



In this section we demonstrate the numerical simulation by the one-dimensional model 

to confirm our analytical expression. This model enables us to simulate with high accuracy, 
since the numerical errors of the algorithm include only round-off errors. Here we investi- 
gate the time evolution of a plane-symmetric density perturbation in an Einstein-de Sitter 
universe model by using a semi-analytical method. In the one-dimensional system many 
plane-parallel sheets move only in perpendicular direction to the surface of these sheets. 
When two sheets cross, they are allowed to pass through each other freely. In this sheet 
system, there is an exact solution until two sheets cross over as follows (Zel'dovich 1970; 
Sunyaev & Zel'dovich 1972; Doroshkevich et al. 1973; Zentsova & Chernin 1980; Yano & 
Gouda 1998): 



(18) 
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Numerical test in the one-dimensional model 



X 



V 



q + B,{t)Si{q) + B2{t)S2{q), 
B,{t)Si{q) + B2{t)S2{q) , 



(19) 
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where q and x are the Lagrangian and the (comoving) Eulerian coordinates, respectively. 
Here, 5*1 (g), and 5*2 (g) are arbitrary functions of q. Bi{t), and B2{t) are the growing mode 
and the decaying mode of hnear perturbation solutions, respectively. We are considering 
the Einstein-de Sitter background universe model, i.e., Bi(t) = a and B2{t) = a~2^ where 
a is the scale factor. The velocity v is the peculiar-velocity normalized by the scale factor. 
We can compute the crossing time of all neighboring pairs of sheets. We use the shortest of 
these crossing times as a time step. Then, we can compute the new positions and velocities 
for all sheets at this crossing time. After two sheets cross, we exchange the velocities of 
the two sheets that just crossed. Then we obtain again Si{q), S2{q), and therefore the exact 
solutions as follows: 

3 2 
Si{q) = -a'^{x - q) + ^a^^w, 

2 2 

S2{q)^-J{x-q)--a-^Jv . (20) 
5 5 

These new exact solutions can be used until again two sheets cross. In this way we obtain 
the exact loci of the sheets by coupling these solutions. This semi-analytic method has good 
accuracy, because we connect exact solutions. 

In our calculation we are going to use 2^^ sheets. Periodic boundary conditions are 
fixed at a length of 2. Here we consider the time evolution from a random-Gaussian initial 
condition with the index of the power spectrum given by 1 or 2. The shape of the triangle 
is chosen as r = 2s. In addition, we take the ensemble average over 100 samples. Figure 1 
shows the power spectrum in the nonlinear gravitational regime where we have chosen the 
initial power index 1. It can be fitted by P{k) oc /c^*^'^^, which agrees with Yano & Gouda 
(1998). From the Fourier transform, the 2PCF obeys ^{r) oc r~°'^^. Figure 2 shows the 
bispectrum in the nonlinear regime when the power spectrum is expressed as in Figure 1. It 
is obvious that the bispectrum is well fitted by B{k,2k) oc k^^'^. Similarly, Figure 4 shows 
the bispectrum in the nonlinear regime when the power spectrum is expressed as in Figure 
3. It is obvious that the bispectrum is well fitted by B{k, 2k) oc These value agree 

with B{k, 2k) oc A;^"*, which is predicted by expression (17) with the relation e = —m/2. In 
other words, these results agree with the hierarchical assumption ^ oc 



5. Summciry 

In this paper we have studied the connection between the 2PCF and the 3PCF in the 
non-linear gravitational clustering regime. Under the scaling hypothesis to the bispectrum 
as well as the power spectrum, we have derived a new formula regarding the 3PCF. It is 
striking that the origin of the power index of the 3PCF has been revealed. Moreover we 
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have checked the vahdity of the analytical expression using numerical simulations in the 
one-dimensional Einstein-de Sitter model. 

Under the scaling hypothesis, we have obtained the connection (17), in which there 
remain two characteristic possibilities for evolution of the phases of the Fourier modes: 
they remain correlated or tend to become uncorrelated in the strongly non-linear regime. 
Expression (17) states that the power index of the 3PCF generally consists of information 
from the third-order statistics, e, as well as that from the 2PCF. 

It is generally believed that the power index of the power spectrum in the non-linear 
regime, m, is decided by the initial power spectrum. Combining this fact and our formula, 
Eq. (17), we can claim that the power index of the 3PCF originates from the initial power 
spectrum. Therefore, it will be important future work to clarify the relationship between the 
power spectrum in the non-linear regime and the initial one. Furthermore, the confirmation 
of the scaling of the 3PCF in three-dimensional N-body simulations is very important future 
work, although it is much more difficult to calculate the 3PCF in three-dimensional space 
than in one-dimensional. 

The radical progress of recent observations needs increasingly advanced and practical 
formalism for analyses containing the higher-order statistics. Hence our new formula is 
expected to be a good indicator to understand the non-linear gravitational clustering. For 
this purpose it will be a subsequent task to check if the new formula is consistent with 
the BBGKY hierarchy. In addition, it will also be important to complete the formula by 
investigating the coefficients as well as the power indices in order to inspect the adjustment 
with the "hierarchical form" and to understand the physical origin of Q. We hope the new 
formula becomes a very useful instrument for analysis of large-scale structure. 

We would like to thank Kei-ichi Maeda, Masahiro Morikawa, and Arihiro Mizutani for 
useful comments. And we are grateful to Takahiko Matsubara and Tsutomu T. Takeuchi for 
many valuable comments and discussions. H.K. is supported by JSPS Fellowship for Young 
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A. Bispectrum under the scaling hypothesis 

Firstly we decompose the absolute value of the Fourier spectrum of density fluctuation 
and the Fourier phase into their averages and the discrepancies: 

Xa^{Xa) + 5Xa^\5tk\, (Al) 

X^ = {X^) + 5Xb=\5tk'l (A2) 
Xc = {Xc) + 5Xc = |5_tfc-tfe'|, (A3) 

Then the bispectrum is rewritten by 

B{tk, tk') = {XaX,X,Y) = {{{Xa) + 5Xa){{Xk) + 5X,){{X,) + 5X,){{Y) + 5Y)) 

= {X,) {X,) (y ) + {X,) {X,) {5Y) + • • • 
{5XM) + --- 
+{Xa) {5X^5X,5Y) + • • • 

+ {5Xa5X^5X^5Y) . (A5) 
Next we rewrite these averages using the variances, Ua, CTft, and Gy and the correlation 

coefficients, Tab, ray, Vabc, Taby and Vabcy, 



2 2 

(7y, 



{Xa) 

{5Y^) 

(SXJXb) = s'^aaCrbVab, 

{SXaSY) = SaaCTyVay, 

{6XaSXb6Xc) = S^GaObOcTabc-, 

{SXaSXf,SY) = s'^aaabayTaby, 
{6XaSXbSXcSY) = s^aaabdcdyrabcy, (A6) 

where s is a constant of order 1. Then 

{XaXhXcY) = aaOb(Tc{Y) + saa(Jb(Jc(Jy{ray + riyy + r^y) + s^{Y)aa(Jb(Jc{rab + + r^a) 

+S^aa(7b(7cay{raby + ^bcj/ + ^cay) + S^(i^)(7a(76Crc?^a6c + aaabacayTabcyi^^) 



This preprint was prepared with the A AS macros v5.2. 
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Now we assume that these constants obey power laws as follows: 

Ua oc (F) oc r'^, (A8) 

Uy^r'^, rab(Xt-"'\ ray(Xr^^\ rabc^t-'^\ raby(Xt-'^^\ rahcy^t'^*'- (A9) 

Substituting (A8) and (A9) into (A7), we find 

{XaXi,X^) (X ^3m/2-a _|_ ^3m/2-/3-72;, _|_ ^3m/2-a-72 _|_ ^3m/2-/3-73j, _|_ ^3m/2-a-73 _|_ ^3m/2-/3-74j^ 

(AlO) 

In order for the 3PCF not to diverge, these scaling factors, a, /3, and 7 should be non-negative 
values. Then the bispectrum can be approximated as 

{XaX^X.Y) oc (All) 

where e is defined as 

e = min{a, (5 + ^2y, /5 + 732/, /5 + 7%)- (A12) 
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Fig. 1. — Power spectrum P{k). We have chosen the initial power spectrum as the power 
law with index n = 1. The power index of the power spectrum in the non-linear regime is 
fitted to -0.75. 
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Fig. 2. — The bispectrum in the nonhnear regime when the power spectrum is expressed as 
in Fig.l. The dashed hne shows a power law with a power index of —1.5. 
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Fig. 3. — The same as Figure 1, but here the initial power index is 2. The power index of 
the power spectrum in the non-hnear regime is fitted by —0.65. 



-16- 




0.5 1 1.5 2 2.5 3 

Log [k] 

Fig. 4. — The bispectrum in the nonhnear regime when the power spectrum is expressed as 
in Figure 3. The dashed hne shows a power law with a power index of —1.3. 



